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1.  Introduction 


Turbulence  is  a  small-scale  phenomena  which  cannot  be  resolved  in  large  or  mesoscale 
weather  forecast  models,  but  which  affects  the  overall  forecasts.  The  problem  that  faces 
the  modeler  is  how  to  parameterize  these  turbulent  motions  within  the  constraints  of  time, 
money  and  computer  power,  while  at  the  same  time  resolving  scales  many  magnitudes 
larger. 

The  modelling  of  turbulence  in  the  boundary  layer  and  in  the  atmosphere  above 
remains  one  of  the  more  difficult  problems  in  the  atmospheric  sciences.  Sometimes 
turbulence  and  diffusion  have  been  construed  as  the  same  process.  In  reality  the  processes 
involve  different  time  and  space  scales  and  are  physically  different.  Consequently,  the 
first  attempts  to  model  turbulence  (Boussinesq,  1877)  took  the  diffusive  formulation  for  the 
molecular  viscosity  of  air,  enhanced  it,  and  called  it  an  eddy  viscosity.  By  specifying  the 
eddy-diffusion  coefficient  in  terms  of  known  variables,  the  hydrodynamical  equations  of 
motion  can  be  closed  and  solved.  Attempts  to  find  suitable  parameterizations  for  the  eddy 
viscosity  based  on  bulk  or  large  scale  variables  has  occupied  a  predominate  role  in  early 
meteorological  boundary  layer  research.  Following  upon  Boussinesq's  constant  eddy 
viscosity  and  Prandtl's  (1925)  mixing  length  theory,  many  complicated  turbulence  models 
have  been  proposed.  In  addition  to  K  theory's  first  order  closure  scheme  (Louis,  1979) 
there  are  higher-order  closure  formulations  (Zeman,  1981;  Wyngaand,  1982;  Mellor  and 
Yamada,  1982;  Andre,  et  al.,  1987)  and  spectral  theories  (Heisenberg,  1948)  for  turbulence. 

In  this  report  we  turn  our  attention  to  the  newly  developed  transilient  turbulence 
parameterization  (Stull,  1984).  The  nonlocal  transilient  turbulence  parameterization  has 
been  compared  with  the  above  mentioned  local  approaches  by  Stull  0984, 1986)  and  to 
turbulent  adjustment  procedures  by  Stull  and  Hasegawa  Q984).  In  many  one  dimensional 
tests  (Stull  and  Driedonks,  1987;  Stull  and  Kraus,  1987)  this  new  formulation  has  yielded 
outstanding  results. 

The  primary  purpose  of  our  investigation  is  to  study  the  importance  and  predictablitily 
of  turbulence  in  a  three  dimensional  numerical  model.  In  our  testing  we  utilized  the 
regional  or  limited-area  Penn  State/NCAR  mesoscale  model.  Our  philosophy  is  to  replace 
separate  formulations  for  K  theory  ,  dry  convective  adjustment,  boundary-layer  evolution 
and  cumulus  parameterization  with  one  formulation,  viz.,  the  transilient  turbulence.  To 
accomplish  our  goals  the  Penn  State/NCAR  regional  model  had  to  be  reformulated  with  the 
transilient  turbulence  parameterization  scheme.  All  existing  K  theory  coding  was 
dismantled  along  with  the  dry  and  moist  convective  adjustment  routines  and  the  boundary- 
layer  parameterization.  A  new  surface  flux  procedure  that  was  compatible  with  the  new 
transilient  turbulence  parameterization  was  installed.  It  allows  for  a  thin  molecular  layer 
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near  the  slab  surface. 

We  selected  the  version  of  the  Penn  State/NCAR  three-dimensional  hydrostatic 
primitive  equation  mesoscale  prediction  model  having  fifteen  vertical  sigma  levels  (Anthes 
and  Warner,  1978).  In  this  model  the  horizontal  grid  included  61  x  46  points  on  a  Lambert 
conformal  mapping  using  a  horizontal  grid  spacing  of  Ax=Ay=80  km.  Other  model 
configurations  include  a  nudging  horizontal  boundary  condition  and  a  two  layer  surface 
slab  formulation.  To  remove  numerical  noise  a  6th  order  implicit  tangent  filter  (Raymond, 
1988)  was  installed.  This  filter  is  very  selective  and  enables  numerical  noise  between  2Ax 
and  4 Ax  to  be  removed  without  altering  meaningful  meteorological  scales.  This  allows  the 
modeler  a  greater  degree  of  freedom  in  testing  turbulence  parameterization  schemes.  The 
time  step  used  is  120s. 

Several  72  h  forecasts  have  been  ran  and  verification  made  for  the  OSCAR  IV  case 
(0000  UTC  22  April  to  0000  UTC  25  April,  198 1 ,  prepared  by  NCAR,  see  Errico  and 
Baumhefner,  1987).  This  data  set  represents  a  spring-time  frontal  situation  in  which 
cyclonic  activity  is  located  initially  over  the  Dakotas  and  an  upper- level  trough  exists  just  to 
its  west  This  system  propagates  eastward  and  intensifies,  but  several  small  scale  waves 
are  also  present  (Errico  and  Baumhefner,  1987).  Additional  cases  will  be  tested  in  the  near 
future  but  only  the  OSCAR  case  is  presented  in  this  first  year  report. 

In  the  next  section  (2)  we  highlight  the  transilient  turbulence  parameterization.  In 
section  3  we  describe  a  surface  flux  scheme  which  contains  a  molecular  layer  modification 
while  in  section  4  the  nature  of  the  6th  order  implicit  filter  is  examined.  Results  from  our 
calculations  arc  presented  in  section  5.  Section  6  summarizes  our  findings  to  date  and  talks 
about  future  activities. 


2.  Transilient  Turbulence 

a.  The  nonlocal  nature  of  transilient  turbulence 

The  transilient  turbulence  parameterization  is  a  nonlocal  first-order  closure  method  that 
models  the  subgrid-scale  mixing  in  the  vertical  (see  Fig.  1).  By  nonlocal,  we  mean  that 
mixing  can  occur  between  all  pairs  of  grid  points  in  a  column  during  any  one  timestep, 
even  if  the  points  are  not  immediate  vertical  neighbors.  Into  each  grid  box  of  interest  i,  the 
fraction  of  air  mixed  from  any  other  grid  box  j,  during  timestep  At  is  represented  by  the 
transilient  mixing  coefficient  Cjj(t,At).  The  amount  of  tracer  S  carried  with  the  air  from 
box  j  to  i  is  the  product  Cjj(t,At)  Sj(t). 
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Fig.  1.  Schematic  showing  grid  boxes  and  flux 
levels.  Only  those  eddies  of  size  3Az  that 
contribute  to  the  flux  across  level  4  are  shown. 

The  set  of  all  coefficients  describing  mixing  between  every  pair  of  points  during  one 
timestep  is  a  matrix  called  the  transilient  matrix.  The  elements  in  each  row  and  each 
column  of  the  transilient  matrix  must  sum  to  one,  to  insure  air-mass  and  tracer-mass 
conservation.  No  element  can  be  negative  without  violating  the  increase  in  entropy  of 
mixing,  thus  each  element  must  satisfy  0  ^  cy  £1. 

The  forecast  equation  for  a  conservative  tracer,  considering  only  turbulence,  is  thus 
described  by  matrix  multiplication  (Stull  and  Driedonks,  1987): 


ii 

Sj(t-t-At)  =  Cjj(t,At)  Sj(t). 


This  same  equation  is  used  to  model  turbulent  mixing  of  heat,  moisture,  and  momentum. 

Fig.  2  shows  how  transilient  matrices  can  be  interpreted:  the  column  index  indicates 
from  where  the  air  is  coming,  and  the  row  index  indicates  to  where  it  is  going.  The 
element  highlighted  in  Fig.  2a  is  in  column  4  row  2,  and  represents  the  mixing  from  box  4 
to  box  2  shown  in  Fig.  2b.  The  arrow  drawn  in  Fig.  2b  is  not  a  physical  eddy,  but  instead 
represents  the  mixing  process  associated  with  a  variety  of  eddies. 
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Fig.  2.  Schematic  (a)  highlighting  one  element  of  a 
transient  matrix,  and  (b)  indicating  the 
associated  mixing  process  between  grid 
boxes  in  a  column. 

b.  The  closure  procedure 

Stull  and  Driedonks  0987)  presented  a  closure  procedure  based  on  the  turbulent  kinetic 
energy  equation.  In  this  nonlocal  closure  scheme  let  Yjj  denote  the  'mixing  potential' 
between  two  boxes  i  and  j.  Yjj  is  the  potential  for  no  mixing  outside  of  grid  box  i.  Then 
the  transilient  coefficients  Cjj  are  defined  by 


Cjj  =  Yjj/  II Y  II  for  i  *  j  , 


with 


cii  =  1  *  X  Cjj 

j*» 


(2a) 


(2b) 


Here  n  denotes  the  number  of  grid  points  in  the  vertical  (15  in  the  Penn  State/NCAR 
model) 

Following  Stull  and  Driedonks  (1987)  the  mixing  potential  is  defined  from  the  turbulent 
kinetic  energy  equation  using  three  scaling  parameters;  a  timescale  of  turbulence  T0,  a 
dimensionless  parameter  1^.  relating  buoyancy  to  shear,  and  a  dimensionless  dissipation 
factor  D,  giving, 

(3) 

Yjj  =  ToAt/(A2)ij2[(AU)ij2  +  (AV)jj2  -  (g/Rc0vi)  (A0v)jj  (Az)jj))  -  DAt/T0  , 


for  i  #  j.  Here  U  and  V  are  the  velocity  components,  0V  the  virtual  potential  temperature, 
g  the  gravitational  constant  while  Azjj  is  the  distance  between  the  centers  of  boxes  i  and  j. 

In  our  calculations  T0= 1000s,  R^.21  and  D=l.  These  values  follow  the  suggestions 
given  in  Stull  and  Driedonks  0987).  They  found  that  T0  should  be  between  100  and  1000 
s.  The  lower  value  was  best  for  neutral  conditions  in  the  boundary  layer.  Rc  represents  a 
critical  Richardson  number  above  which  no  turbulence  takes  place.  D  makes  allowances 
for  molecular  viscosity  mixing. 

The  internal  mixing  potential  for  box  i  defined  by  Yfi  is  restricted  to  be  as  large  as  any 
of  the  elements  in  the  same  row.  To  avoid  convective  overturning  we  furthermore  restrict 
that  the  mixing  potential  be  nonincreasing,  i.e.,  away  from  the  main  diagonal  the  mixing 
potential  does  not  increase  in  magnitude  within  any  row  of  the  transient  matrix. 

c.  The  numerical  procedure 

During  each  timestep,  the  above  1-D  transilient  method  is  applied  to  each  column  of  the 
grid  points  in  the  model.  In  addition  to  performing  the  vertical  mixing,  the  locations  of 
turbulence  (both  in  the  boundary  layer  and  aloft  as  clear  air  turbulence)  are  saved.  These 
locations  are  then  output  as  CAT  regions.  In  addition,  horizontal  mixing  (proportional  to 
the  depth  of  the  vertical  mixing)  is  applied  to  only  those  neighboring  grid  points  that  were 
turbulent.  The  net  result  is  a  qilasi-3-D  turbulence  scheme  that  automatically  forecasts 
boundary-layer  evolution  and  patchy  CAT. 

3.  Heat  Flux  With  a  Molecular  Layer  Included 

The  heat  flux  used  in  the  surface  layer  turbulence  calculations  within  the  Penn  State/ 
NCAR  model  is  computed  from 

Hs  =  -pCpku^T*  2  -Cj(0a  -  0S) ,  (4) 

where  T*  is  given  by 

T*  =  (  0  -  Qjnnizjz^  -  yH)  •  (5) 

The  parameter  p  is  the  density,  k  is  the  von  Karman  constant,  Cp  is  the  specific  heat  of  air 


yuuu  vx  w  vx  vy  i/xvvvvyv .' 


at  constant  pressure,  6(  is  the  potential  temperature  at  the  nominal  height  z,  while  0,  is  the 
corresponding  air  potential  temperature  at  the  lowest  model  level  z,.  The  friction  velocity, 
u*,  is  expressed  as 

u*=Max(kV/ln(za/zo  -  V) ,  u*  ) ,  (6) 

M  ° 

while  over  water  the  definition  of  u*  is  changed  to 

u*=(4  x  10-4  V2-55)0J  ,  (7) 

where  u*o  =.l  ms'^  and  V=(Va2  +  Vg2)0-5.  Here  Va  is  the  wind  speed  at  height  za  and 
Vc  is  a  convective  velocity  (Anthes  et  al.,  1987).  The  roughness  height  is  zQ.  The 
nondimensional  stability  parameters  yH  and  \|/M  are  a  function  of  the  bulk  Richardson 
number  (Anthes  et  al.,  1987). 

As  will  be  shown  below  it  is  very  important  to  include  the  molecular  layer  in  the  heat 
flux  calculations.  The  molecular  layer  is  of  depth  z^ .  The  molecular  heat  flux  at  is 
defined  by 

HSn  =  -pKH^(Ts  -  Tg)^  *  -c2(0s  -  0g).  (8) 

Here  Kh^  is  the  molecular  conductivity  in  air. 

Solving  for  0S  in  eq.  (8),  substituting  into  Eq.  (4)  and  requiring  the  molecular  and 
turbulent  heat  fluxes  to  be  identical  in  the  surface  layer,  i.e.,  H  *  Hs  =  Hs(i ,  yields 

H=  -c,(0a  -  0g)(l+  Cl/c2)  i  =  -CiC3(0a  -  0g)  (9) 

The  ratio  c^Cj  is  defined  as 

ci/c2  =  cp  k2  v  Zp /[KHplnfza/zo  -  v^^nCZg/z,,  -  \|/H)]  (10) 

For  typical  values  of  V,  za,  z^,  z0 ,  e.g.,  V=5ms1.  za=50  m,  Zp=.01  m  and  z0  between 
0.01  and  0.5,  we  obtain  for  c3  *  (1+  cj/c^'1  a  range  of  values  between  0.20  and  0.06. 
Consequently  it  is  clear  that  the  molecular  layer  reduces  the  heat  flux,  Eq.  (9),  by  at  least 
80%.  Its  inclusion  into  the  model  shelters  the  calculations  from  unrealistic  gradients 
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associated  with  very  steep  super  adiabatic  layer  found  at  the  surface  during  the  day,  and 
from  the  reverse  situation  at  night  In  our  model  runs  we  have  fixed  the  value  of  C3  at 
0.2. 

4.  The  Implicit  Tangent  Filter 

The  set  of  symmetric  low-pass  implicit  tangent  filters  (Raymond,  1988)  of  order  2p  are 
defined  by 

[S2p]u„f  +  (-1)P£[L2p]uF=  [S2p]Un  .  (11) 

Here  u„F  is  the  filter  variable  while  u^  is  the  unfiltered  quantity.  In  Eq.  (11)  L2u„  is  the 
finite-difference  analog  of  the  second  derivative  of  u  at  grid  location  n  multiplied  by  the 
square  of  the  grid  step  size,  i.e.,  L2^  =  u^i  -  2u„  +  u„+1.  Similarly,  order  2p  of  the 
operator  is  analogous  in  the  same  way  to  a  2p  derivative.  The  coefficients  associated  with 
the  (L)2?  operation  are  identical  to  those  in  the  binomial  expansion  of  (a  -  b)2P.  Also 
S2un=  Un_j  +  2un  +  Uji+i  ;  consequently  (S/2)2  may  be  thought  of  as  an  averaging  operator. 
The  coefficients  for  (S)2p  are  the  2p  row  entries  in  Pascal's  triangle. 

The  implicit  formulation  requires  the  inversion  of  a  banded  diagonally  dominant  matrix 
to  determine  the  filtered  set  u„F  for  all  n.  For  p=l  the  filter  utilized  in  Pepper,  et  al.,  (1979) 
is  recovered.  The  filter  parameter  e  is  varied  to  control  the  degree  of  filtering,  i.e.  E=0 
gives  no  filtering  while  E=1  will  reduce  a  4Ax  feature  by  fifty  percent  for  any  order  2p. 

The  characteristics  of  the  low-pass  implicit  tangent  filter  are  best  illustrated  by  examining 
the  amplitude  response  function. 

Assuming  a  Fourier  expansion  for  both  Uj,  and  u„F  in  Eq.  (1 1)  and  forming  the  ratio  of 
the  amplitudes  yield,  after  simplification,  an  amplitude  response 

F(k)  =  [  1  +  e  tan2P(K8/2)]-] .  (12) 

Here  5  is  the  grid  step  distance  while  K  is  the  wave  number.  This  response  is  the  same  as 
that  given  by  the  recursive  tangent  filters  described  in  Otnes  and  Enochson  (1978). 

However  the  implicit  tangent  filter  is  nevertheless  much  easier  to  work  with  because  the 
coefficient  weights  are  known  and  boundary  conditions  encountered  in  limited  area 
modelling  are  easily  handled  in  the  implicit  formulation.  Tedious  calculations  are 
otherwise  required  to  determine  the  coefficients  of  the  recursive  tangent  filter  (Otnes  and 
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Enochson,  1978)  and  boundary  conditions  are  not  obvious. 

In  Fig.  (3)  the  response  of  the  filter  is  shown  after  24  hrs  (720  applications  with 
£=.0075)  and  is  compared  with  the  4th  order  K  theory  horizontal  smoothing  normally  used 
in  the  Penn  State/NCAR  regional  model.  Note  for  wave  number  6Ax  that  the  responses 
are  much  closer  to  the  unfiltered  value  of  1  when  using  the  6th  order  implicit  filter,  whereas 
the  old  smoothing  operator  removes  almost  all  of  these  feature,  especially  when  the 
maximum  K  value  is  used. 
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Fig.  3.  The  filter  response  verses  wave  length  nAx 
after  720  applications  for  the  implicit 
tangent  filter  and  for  the  Penn-State/NCAR 
4th  order  K  theory  diffusion. 

In  the  model  forecasts  the  filter,  with  £=.0075,  is  applied  at  every  time  step  to  the 
horizontal  wind  velocity  components,  the  temperature  and  the  mixing  ratio  fields.  The  filter 
parameter  e  was  selected  to  give  the  lowest  acceptable  smoothing  while  still  maintaining 
stability  and  reasonable  smooth  fields.  The  application  of  the  filter  takes  place  after 
completion  of  the  new  time  step  and  following  the  turbulence  parameterization. 
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5.  Discussion  of  Results 

In  Fig.  4  we  illustrate  for  the  OSCAR  IV  case  the  percent  of  the  total  number  of 
horizontal  grid  points  that  contain  turbulent  exchanges,  as  a  function  of  pressure  and  time, 
as  predicted  from  a  72  hr  forecast  with  the  transilient  turbulence  parameterization.  In  Fig.  4 
note  that  above  the  500  mb  level  approximately  one  percent  or  less  of  the  total  number  of 
the  horizontal  grid  points  are  turbulent  The  maximum  at  high  levels  occurs  during  the 
strongest  cyclonic  intensification  which  takes  place  between  hours  18  through  48.  In  the 
lowest  levels  of  the  boundary  layer  more  than  75%  of  the  total  horizontal  grid  points  are 
turbulent  during  the  peak  heating  in  the  diurnal  cycle.  Clearly  there  is  a  significant  diurnal 
variation  in  the  boundary  layer. 
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Fig.  4.  The  percent  of  the  grid  points  turbulent  as  a 
function  of  time  and  pressure. 

The  horizontal  distribution  of  turbulence  in  the  boundary  layer  is  shown  in  Figs.  5a,b. 
In  Fig.  5a  the  distribution  of  the  turbulence  21  hours  into  the  forecast  (near  3  pm  local  time) 
at  0=0.94  corresponds  closely  with  the  heated  land  mass.  Only  in  the  midwest  (associated 
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with  the  major  cyclonic  activity)  is  there  a  turbulence  free  zone.  This  most  likely  occurs 
because  of  die  stabilizing  radiative  and  evaporative  properties  associated  with  the 
precipitating  clouds.  At  higher  levels  in  die  atmosphere  (c=.74)  the  turbulence  is  confined 
to  the  mountainous  regions  as  shown  in  Fig.'Sb . 
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Fig.  5.  The  horizontal  distribution  of  turbulence  at  (a) 

o=.94  and  (b)  at  a=.74  twenty  one  hours  into  the 
forecast. 
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Vertical  cross  sections  through  the  center  of  our  region  (taken  west  to  east)  are 
presented  in  Figs.  6a,b  showing  die  potential  temperature  0  (K)  and  the  mixing  ratio 
(g  kg'1).  Regions  with  turbulence  are  within  or  under  the  wide  solid  line  found  in  the  0 
field  in  Fig.  6a.  The  turbulent  boundary  layer  appears  deepest  over  the  mountains  and  just 
ahead  of  the  front  as  shown  between  grid  points  30  and  45,  which  is  very  realistic.  Some 
turbulence  is  also  occurring  at  mid-levels  within  the  frontal  zone  above  grid  points  29  and 
30.  Note  that  some  folding  in  the  0  contours  occur  at  mid  levels  between  grid  points  25 
and  35  indicating  the  presence  of  the  cold  front.  Clearly  the  largest  amount  of  moisture  is 
found  in  the  warm  sector.  This  is  reflected  in  the  contours  of  mixing  ratio  in  Fig.  6b 
which  have  their  greatest  vertical  extent  between  grid  points  30  and  45. 


0  10  20  30  40  SO  60 

Grid  Point 


(b) 


0  10  20  30  40  50  60 

Gfid  Point 


Fig.  6.  Vertical  cross  sections  of  the  (a)  potential 

temperature  (K)  and  (b)  mixing  ratio  (g  kg'1). 
Turbulence  indicated  by  the  thick  line  and 
stipple  in  (a). 


Fig.  7  shows  the  rms  error  at  850  mb  for  the  forecast  temperature.  Results  from  four 
forecasts  are  displayed.  The  control  run  is  the  Penn  State/NCAR  model  with  the  existing 
Blackadar  boundary  layer  scheme  and  K  theory  horizontal  diffusion.  The  control  was  ran 
with  a  Kuo  cumulus  parameterization  and  with  an  explicit  cloud  scheme,  labeled  (clouds)  in 
Fig.  7.  The  transilient  turbulence  parameterization  version  was  ran  similarly  except  no 
cumulus  parameterization  was  utilized  in  any  of  these  computer  simulations.  Note  in  Fig. 

7  that  the  rms  errors,  determined  when  forecast  and  radiosonde  values  are  compared,  are 
least  for  the  transilient  approach  with  explicit  clouds,  except  near  the  end  of  the  72  hr 
forecast  Overall,  carrying  the  clouds  explicitly  made  little  difference  in  these  rms  statistics. 


RMS  error  in  850  mb  forecast  temperature  .vs.  radiosonde  reports 
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Fig.  7. 


The  rms  error  in  the  850  mb  temperature  for  the 
control  forecast,  for  the  standard  Penn  State 
/NCAR  package  of  the  Blackadar  boundary  layer 
formulation  with  K  theory,  with  the  Kuo 
cumulus  parameterization  and  with  explicit 
clouds.  Also  illustrated  is  the  transilient 
results  with  and  without  explicit  clouds. 
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The  affect  of  explicit  clouds  is  however  very  clear  in  the  mean  error  in  the  850  mb 
temperatures  shown  in  Fig.  8.  Note  that  a  change  in  sign  is  associated  with  whether  the 
clouds  are  carried  explicitly  are  not  The  transilient  turbulence  scheme  has  zero  mean 
error  at  hours  12  and  24  ,  and  has  the  smallest  mean  error  of  the  four  forecasts.  The 
control  run,  with  and  without  explicit  clouds,  generally  has  much  larger  mean  errors. 


Mean  error  in  850  mb  forecast  temperature  .vs.  radiosonde  reports 
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Fig.  8.  Mean  error  for  the  cases  displayed  in  Fig.  7. 


Evidence  is  presented  in  Fig.  9  which  shows  that  a  systematic  bias  exists  in  the  mean 
error  for  the  850  mb  geopotential  heights.  Note  that  the  control  run  with  the  Kuo  cumulu: 
parameterization  has  almost  no  mean  error  while  the  same  control  forecast  made  using  an 
explicit  cloud  formulation  has  the  largest  error.  Part  of  this  error  may  be  ascribed  to 
inaccuracies  in  the  forecasted  phase  velocity.  However,  evaporative  and  radiative  feedback 
processes  in  the  surface  calculations  are  also  of  major  importance.  The  error  in  the 


transilient  turbulence  approach  suggests  that  the  low-level  forecast  temperatures  are  too 
warm  (approximately  one  degree).  This  is  confirmed  by  comparing  the  temperatures 
predicted  at  the  lowest  sigma  level  (not  shown).  The  geopotential  height,  being  an  integral 
of  the  virtual  temperature  retains  errors  that  originate  at  the  lower  levels.  In  this  regard  we 
must  scrutinize  our  new  surface  flux  calculations  and  the  surface  flux  balances  to  make  sure 
we  are  getting  the  best  solution,  otherwise  we  can  introduce  biases  into  our  forecast 
calculations.  Additional  testing  still  needs  to  be  done  in  this  regard. 

Mean  error  in  forecast  850  mb  heights  .vs.  radiosonde  reports 
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Fig.  9.  Same  as  Fig.  8  except  for  the  mean  error  in  the 
850  geopotential  heights  (m). 


In  Fig.  10  the  mean  sea-level  pressure  48  hrs  into  the  forecast  is  shown  for  the  control 
runs,  Kuo  cumulus  parameterization  (10a),  explicit  cloud  (10b),  and  for  the  transilient 
version  with  explicit  cloud  (10c)  and  for  the  verification  analysis  (lOd).  Note  that  the  low 
pressure  centered  over  the  great  lakes  is  best  represented  in  the  transilient  forecast  which 
has  the  correct  pressure  of  997  mbs  but  over  a  reduced  area  as  compared  to  the  verifying 
analysis  in  Fig.  lOd.  The  contraction  of  the  surface  low  in  the  transilient  case  is  similar  to 
the  control  case,  with  explicit  cloud  calculations,  which  has  a  central  pressure  of  996  mb. 
The  control  with  the  Kuo  scheme  has  a  990  central  pressure,  so  the  low  is  deepened  too 
much.  As  a  consequence  of  this  the  1008  pressure  contour  is  however  in  a  location  closer 


to  that  indicated  in  the  verifying  analysis.  Otherwise  over  most  of  the  remaining  area  the 
control  cases  are  slightly  better  by  a  small  margin,  but  the  presence  of  topography  is  a 
complicating  factor  making  interpretation  difficult. 
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Fig.  10.  Mean  sea  level  pressure  for  (a)  the  control  run  with  Kuo 

cumulus  parameterization  and  (b)  with  explicit  clouds.  The 
transilient  turbulence  parameterization  scheme  with  exolidt 
clouds  is  shown  in  (c).  The  verifying  48  h  analysis  is 
given  in  (d).  Contour  interval  is  4  mb. 
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Fig.  10.  Mean  sea  level  pressure  for  (a)  the  control  run 
with  Kuo  cumulus  parameterization  and  (b)  with 
explicit  clouds.  The  transilient  turbulence 
parameterization  scheme  with  explicit  clouds  is 
shown  in  (c).  The  verifying  48  h  analysis  is 
given  in  (d).  Contour  interval  is  4  mb. 
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To  gain  some  idea  of  the  response  at  the  surface  we  now  turn  our  attention  to  Figs.  1  la 
and  lib  showing  the  surface  slab  or  skin  temperature,  48  hrs  into  the  forecast  at  0000  UTC 
24  April,  for  the  control  case  with  Kuo  cumulus  parameterization  and  the  equivalent 
transilient  turbulence  simulation.  A  comparison  of  these  two  figures  shows  that  they  are 
nearly  the  same  in  the  north-eastern  portion  of  our  region  while  in  the  western  part  the  skin 
temperature  differences  are  up  to  4  degrees  is  some  locations.  We  expected  and  wanted 
slightly  warmer  surface  skin  temperatures  in  the  transilient  version.  Remembering  that  the 
fluxes  are  being  reduced  by  C3  in  Eq.  (9),  which  compensates  for  the  molecular  layer, 
means  that  the  overall  heat  flux  felt  by  the  lowest  layer  should  be  about  the  same  as  that 
observed  in  the  control  case. 
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Fig.  11.  The  forecasted  skin  temperature,  at  hour  48,  for 
the  control  (a)  with  the  Kuo  cumulus  scheme, 
and  for  (b)  the  transilient  turbulence  formulation. 
Contour  interval  is  3  K. 
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Rainfall  from  the  control  and  transilient  72  hr  forecasts  have  about  the  same  bias  in  all 
categories,  as  shown  in  Fig.  12,  except  for  the  last  (2.64  cm)  or  largest  rainfall 
accumulation.  In  the  latter  case  note  that  all  forecasts  had  a  large  positive  bias.  We 
believe  this  bias  was  enhanced  for  the  transilient  turbulence  parameterization  version 
because  of  two  factors.  The  removal  of  K-theory  horizontal  diffusion  increases  the  rainfall 
by  the  elimination  of  entrainment,  because  there  is  no  horizontal  mixing  except  where  there 
is  vertical  mixing.  Fig.  5  showed  a  significant  turbulence-free  zone  in  the  boundary  layer, 
induced,  we  believe,  because  of  evaporative  and  radiative  cooling  associated  with 
precipitating  clouds  (based  on  precipitation  plots  not  shown).  Additionally,  the  transilient 
approach  tends  to  enhance  boundary  layer  winds.  Consequently,  additional  moisture  is 
transported  into  the  warm  sector  by  the  southerly  component  of  the  horizontal  wind  (not 
shown). 

Bias  scores  (72  hrs)  with  and  without  explicit  clouds 
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Fig.  12.  Bias  in  the  rainfall  (cm)  for  the  control,  with  Kuo 
or  explicit  clouds  and  for  the  transilient 
turbulence  with  or  without  explicit  clouds. 

Occasionally  models  have  trouble  with  excessive  condensation  and/or  overdevelopment 
as  discussed  by  Zhang,  et  al.  (1988).  Such  over  development  results  in  and  is  assisted  by 
positive  feedback  mechanisms  in  latent  heat  release,  moisture  convergence  and  surface 
pressure  decreases.  Including  realistic  physics  like  evaporation  and  water  loading  reduces 
this  problem  (Zhang,  et  al.,  1988).  In  our  examination  of  the  mean  sea  level  pressure, 

Fig.  10,  we  saw  that  the  transilient  turbulence  parameterization  scheme  did  not  cause 
excessive  development  in  low  pressure.  Including  the  water  loading,  via  explicit  clouds, 
however  did  increase  the  accuracy  of  predicting  the  mean  sea  level  pressure  for  both  the 
transilient  parameterization  approach  and  for  the  control  forecast. 
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6.  Conclusions  and  Future  Activities 


We  are  very  pleased  that  the  transilient  turbulence  parameterization  scheme  does  not 
cause  the  Penn  State/NCAR  regional  model  to  become  numerically  unstable.  We  are  also 
pleased  with  some  of  the  eariy  findings.  In  general  we  find  that  the  new  calculations  are 
very  similar  to  the  control  simulations.  We  believe  that  our  statistics  will  improve  when 
we  have  the  new  surface-layer  flux  scheme  tuned  correctly.  Currently  the  lowest  level 
temperature  values  are  not  now  as  accurate  as  they  should  be.  This  introduces  a  bias  into 
the  geopotential  heights,  as  discussed  in  section  5,  which  gets  carried  throughout  the  entire 
depth  of  the  model  because  the  heights  are  calculated  by  an  integral  process.  Our  testing 
has  taken  somewhat  longer  than  anticipated  because  of  the  number  of  parameters  in  the 
existing  Penn  State/NCAR  surface  flux  parameterization  to  be  tuned.  We  plan  to  complete 
all  tests  within  the  next  two  months. 

We  plan  to  examine  the  differences  between  standard-model  forecasts  made  with 
cumulus  parameterization  and  explicit  cloud  calculations  against  our  improved  model  with 
and  without  explicit  clouds.  A  detailed  study  of  this  type  should  be  quite  relevant  and 
important  for  future  attempts  to  adapt  the  transilient  turbulence  technique  to  cloud 
parameterizadon. 

To  date  out  testing  of  the  transilient  turbulence  parameterization  has  utilized  just  one 
data  set,  i.e.,  the  OSCAR  IV  case.  As  soon  as  we  are  satisfied  that  we  are  doing  the  best 
we  can  we  will  begin  simulating  other  cases.  We  plan  to  make  forecasts  for  at  least  six 
independent  cases  so  that  we  can  make  some  definitive  statement  regarding  the  accuracy, 
impact  and  desirability  of  the  transilient  turbulence  parameterization  scheme.  We  also 
want  to  get  at  least  one  data  set  that  is  less  than  five  years  old  so  that  we  can  make 
comparisons  between  our  forecasts  of  turbulence  and  the  repons  of  CAT  obtained  from 
airplane  pilots. 
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